For many years seismic exploration for oil and gas has involved the use of a source of seismic energy and its reception by an array of seismic detectors, generally referred to as geophones. When used on land, the source of seismic energy can be a high explosive charge electrically detonated in a borehole located at a selected point on a terrain, or another energy source having capacity for delivering a series of impacts or mechanical vibrations to the earth's surface. The acoustic waves generated in the earth by these sources are transmitted back from strata boundaries and reach the surface of the earth at varying intervals of time, depending on the distance and the characteristics of the subsurface traversed. These returning waves are detected by the geophones, which function to transduce such acoustic waves into representative electrical signals. In use an array of geophones is generally laid out along a line to form a series of observation stations within a desired locality, the source injects acoustic signals into the earth, and the detected signals are recorded for later processing using digital computers where the data are generally quantized as digital sample points such that each sample point may be operated on individually. Accordingly, seismic field records are reduced to vertical and horizontal cross sections which approximate subsurface features. The geophone array is then moved along the line to a new position and the process repeated to provide a seismic survey. More recently seismic surveys involve geophones and sources laid out in generally rectangular grids covering an area of interest so as to expand areal coverage and enable construction of three dimensional (3D) views of reflector positions over wide areas.
As oil and gas production from known reservoirs and producing provinces declines, explorationists seek new areas in which to find hydrocarbons. Many of the new areas under exploration contain complex geological structures that are difficult to image with 2D techniques. Accordingly, 3D seismic processing has come into common use for mapping subterranean structures associated with oil and gas accumulations. Geophysics, however, are well aware that a 2D seismic record section or 3D view is not a true reflectivity from the earth, but is instead a transformation of the earth's reflectivity into a plane where each recorded event is located vertically beneath the source/receiver midpoint. Steep dip, curved surfaces, buried foci, faults and other discontinuities in subterranean structure each contribute their unique characteristics to the seismic record and, in complexly faulted and folded areas, make interpretation of the geological layering from the seismic record extremely difficult. Migration is the inverse transformation that carries the plane of recorded events into a true 3D reflectivity of the earth, thereby placing reflections from dipping beds in their correct location and collapsing diffractions.
Of the various available migration methods, wave-equation migration is considered to be superior because it is based on accurate propagation of seismic waves through complex models of the earth. Wave equation computations using numerical techniques have led to a procedure called reverse time migration. By this procedure, the wavefield recorded at the surface is imaged in depth using a model of earth velocities in a numerical solution of the wave-equation. The wavefield at the surface is used as a boundary condition for the numerical computations. Proceeding by inserting the data at the surface of a computational grid for each record time step, starting with the last recorded time sample, and ending with the first, the wavefield migrates to the position from which the reflections originated.
A wave equation model is represented by solutions of the second order partial differential scalar wave equation: ##EQU1## where the .gradient..sup.2 operator stands for ##EQU2## V is earth medium velocity, and .PHI. is the physical quantity measured: particle velocity, acceleration, displacement or pressure.
The scalar wave equation can also be written as a set of coupled first order partial differential equations: ##EQU3## Where P is the acoustic pressure, .orgate. is the particle displacement velocity, .rho. is the density of the medium, and c is the propagation speed of the medium.
A practical procedure for doing reverse time migration (RTM) is disclosed in a publication, Mufti, I. R., et al, "Finite-Difference Depth Migration of Exploration-Scale 3D Seismic Data," Geophysics, Vol. 61. No. 3 (May-June 1996), which is incorporated herein by reference. RTM, which requires enormous computer resources as compared to simpler or less accurate migration algorithms, has recently been applied to 3D seismic data. An improved image results from accuracy (dynamic as well as kinematic) of the finite-difference method over conventional normal-movement and raytracing based seismic imaging methods. In this RTM procedure a finite-difference earth model, which is based on the best estimate of subsurface velocities, is required. This involves dividing the model space into a large number of elementary grid blocks, and assigning a velocity value to each grid.
A common problem with the finite-difference migration method is that simulation of wave propagation in an extensive portion of the earth must be modeled with limited computer resources, i.e., mainly limited central RAM memory in the computer such that artificial boundaries are introduced to limit the computational size of the model. Accordingly, artificial boundaries that act as perfect reflectors are imposed by computer memory limitations, and if not properly handled cause unwanted reflection waveforms. These unwanted reflection waveforms may disrupt the image in complex ways and are not easily removed after the modeling and imaging processes are carried out. In two-dimensional modeling, several approaches have been derived which effectively reduce the amplitudes of these artificial boundary reflections so that they do not cause problems in the region of interest. In 3D simulation, however, the approaches effective in two dimensions are insufficient to reduce boundary reflections to an acceptable level. The basic computing requirements for successful 3D finite-difference acoustic wave equation computations include a large memory and a high-speed processor, however, an urgent need exists for a cost effective and straightforward method to reduce artificial boundary reflections in finite-difference 3D wave equation computations.
Accordingly, an object of the present invention is to eliminate reflections from artificial boundaries in 3D finite-difference acoustic wave equation computations for migration of seismic reflections.
A more specific object is to cancel artificial boundary reflections in computer implemented finite-difference calculations without significantly increasing computer memory capacity and/or computation time.
A still further object of this invention is to produce a computer program which generates high resolution images of seismic wave propagation.